function error = posteriori_elliptic(dofs,V,T,d,c,Int_UU_cell,Int_UV_cell,Int_UW_cell,Int_VV_cell,Int_VW_cell,Int_WW_cell)
% error =
% posteriori_elliptic(dofs,V,T,posT,d,c,Int_UU_cell,Int_UV_cell,Int_UW_cell,Int_VV_cell,Int_VW_cell,Int_WW_cell)
% now the estimator is only fit for elliptic problems
nt = size(T,1);
error = zeros(nt,1);
for k = 1:nt
    %comput the local bending matrix
    Int_UU = Int_UU_cell{d(k)};
    Int_UV = Int_UV_cell{d(k)};
    Int_UW = Int_UW_cell{d(k)};
    Int_VV = Int_VV_cell{d(k)};
    Int_VW = Int_VW_cell{d(k)};
    Int_WW = Int_WW_cell{d(k)};
    V1=V(T(k,1),:);V2=V(T(k,2),:);V3=V(T(k,3),:);
    u = (V2(2)-V3(2))^2+(V2(1)-V3(1))^2;
    v = (V1(2)-V3(2))^2+(V1(1)-V3(1))^2;
    w = 2*((V2(2)-V3(2))*(V1(2)-V3(2))+(V2(1)-V3(1))*(V1(1)-V3(1)));
    area = abs((V2(2)-V3(2))*(V1(1)-V3(1)) - (V1(2)-V3(2))*(V2(1)-V3(1)))/2;
    LocK = (u*u*Int_UU - u*w*(Int_UW+Int_UW') + v*v*Int_VV +...
        w*w*Int_WW - v*w*(Int_VW+Int_VW') + u*v*(Int_UV+Int_UV'))/(16*area^3);   % J = 2*tri_area(K)
    % get dofs of this triangle
    loc_dof = dofs(k,1):dofs(k,2);
    c_loc = c(loc_dof);
    % find longest side's length
    sides = [(V(T(k,2),1)-V(T(k,3),1))^2 + (V(T(k,2),2)-V(T(k,3),2))^2; ...
        (V(T(k,1),1)-V(T(k,3),1))^2 + (V(T(k,1),2)-V(T(k,3),2))^2; ...
        (V(T(k,1),1)-V(T(k,2),1))^2 + (V(T(k,1),2)-V(T(k,2),2))^2];
    longest =  max(sides);
    error(k) = sqrt(abs(c_loc'*LocK*c_loc))*longest;
end
    